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Abstract 

The asymptotic frequency, u>, dependence of the dynamic viscosity 
of neutral hard sphere colloidal suspensions is shown to be of the form 
r]oA{(j)){u!Tp)~ 1 / 2 , where A((f>) has been determined as a function of 
the volume fraction 0, for all concentrations in the fluid range, 770 is 
the solvent viscosity and Tp the Peclet time. For a soft potential it is 
shown that, to leading order in the steepness, the asymptotic behavior 
is the same as that for the hard sphere potential and a condition for the 
cross-over behavior to 1/uTp is given. Our result for the hard sphere 
potential generalizes a result of Cichocki and Felderhof obtained at 
low concentrations o and agrees well with the experiments of van der 
Werff et ala, if the usual Stokes-Einstein diffusion coefficient Dq in 
the Smoluchowski operator is consistently replaced by the short-time 
self diffusion coefficient D s (4>) for non-dilute colloidal suspensions. 



Keywords: viscosity, visco-elasticity, rheology, colloidal suspensions, hard- 
spheres, soft-spheres. 



1 Introduction 



The visco-elastic behavior, i.e. the frequency dependent viscosity rj((f),u), 
of concentrated neutral hard sphere colloidal suspensions has been obtained 
in the benchmark experiments of Van der Werff et alo for volume frac- 
tions 0.44 < (f) < 0.58, where 4> = irna^/G, with n the number density of 
the colloidal particles of diameter a. The frequency dependence was found 
to be qualitatively similar to that obtained theoretically by Cichocki and 
FelderhofH for dilute suspensions from an exact solution of the two particle 
Smoluchowski equation for two Brownian particles without hydrodynamic 
interactions. 

An approximate theory for concentrated colloidal suspensions was devel- 
oped by Verberg et alH which agreed well with the experimental results of 
Van der Werff et al for such suspensions. In particular, the asymptotic be- 
havior of the (complex) viscosity for large frequencies uj was given correctly 
as ~ TioA(4>)(l + i)/ yJujTp, where 770 is the viscosity of the solvent, A(4>) an 
amplitude and rp a characteristic Brownian particle interaction time, the 
Peclet time, defined below. However, the amplitude A(<j)) was at small <f> a 
factor two smaller than the exact value obtained by Cichocki and Felderhof 
at low densities and it was too high when compared with the experiments of 
Van der Werff et al at high densities. This difference in asymptotic behav- 
ior did not affect the good agreement with experiments carried out in the 
reduced form used by Van der Werff et al@ @' El . 

In the theory of Verberg et al, r](4>,uj) was obtained as a sum of two 
terms: a short time - infinite frequency - contribution ^^((j)) on the very 
short Brownian-time scale rp(~ 10 _9 s) where the Brownian particle forgets 
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its initial velocity and a long time contribution, on the very much longer 
Peclet-time scale rp(~ 10 _4 s), involving mode- mode coupling contributions 
associated with two cage-diffusion modes, that describe the diffusion of each 
colloidal particle out of the cage in which it finds itself in a concentrated 
colloidal suspension!! : 

r]((j),Lo) = r/oo^) + r/ mc (0,w) (1.1) 
For large u, the mode-mode coupling contributions rj mc ((p,uj) reduces to: 

VrnciM = U 2 x(<Pf /2 ^={i + i)m + O(-) (1.2) 

where x{4>) 1S the equilibrium radial distribution function g eq (r;<j)) at con- 
tact, i.e. x{4>) = 9eq{f = cr;<f>), where r is the distance between two hard 
spheres of diameter a and Tp = <7 2 /4L>o- Here Dq is the Stokes-Einstein 
diffusion coefficient 

Oo = (1.3) 

where kp is Boltzmann's constant and T the temperature of the colloidal 
suspension. 

For low concentrations (ft — ► 0, x(0) — > 1 and r] mc (<j),u>) reduces to: 

7] mc (4>,uj) = r](<f>,u) - Voo(<P) ^ I p — (1 + i)r/o (1-4) 
while Cichocki and Felderhof obtainS: 

VCF(<t>,u) - ^Oft) ^» ^ F — (1-H)r? (1.5) 

5 ^OOTp 

The different coefficient in eq.(1.4) for the approach to 7700 (0) is due to the 
approximate nature of r] mc ((f), cj)@ . 
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The purpose of this paper is to obtain the exact asymptotic behavior 
of rj{(f),uj) for large u for all <fi studied by van der Werff et al@l, i.e. an 
extension of Cichocki and Felderhof's result to high concentrations, as well 
as its behavior for a soft potential. 

In the next section we will give the basic equations. In section 3 we 
will calculate the asymptotic frequency dependent viscosity, for a soft, but 
very steep potential, starting from the Green-Kubo expression. In section 
4 we will give the result for a hard sphere potential, as the limit of a soft 
potential. We will end with a short discussion on the soft potential result. 

2 Basic equations 

In order to obtain the asymptotic behavior of r/((j), to) for large oj for concen- 
trated suspensions we start from a general Green-Kubo expression for the 
frequency dependent viscosity of a colloidal suspension® : 



Eq.(2.1) gives the linear response of the suspension to an applied shear rate 
j(t) = 7oe _,w * with finite frequency u and vanishing shear rate amplitude 
7o — ► 0. In eq.(2.1), j3 = l/Zc^T, V is the volume of the colloidal suspension, 
while p v (t;(J)) is the stress-stress auto correlation function defined by: 



where the brackets denote a canonical equilibrium ensemble average. The 
microscopic stress tensor Y?l (r N ) is given by: 




(2.1) 




(2.2) 



N 




(2.3) 
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with Ti the position of particle i (i = 1,...,N), r N = t\,...,Tn, Fj = 
— Vi&(r N ) the total force on particle i (Vj = <9/<9r;), $(r N ) the total po- 
tential energy of the colloidal particles and: 

N 

n(r N ;<f>)= J2(^i + ^ F i(r N ))-^ij(r N )-^j (2-4) 

the TV-particle Smoluchowski operatorfl HI [I , the colloidal analogue of the 
Liouville operator for atomic liquidsE, with Dy(r^) the diffusion tensor, 
incorperating hydrodynamic interactions. This diffusion tensor determines 
the velocity imparted to particle i by a force acting on particle j. In the 
absence of hydrodynamic interactions, i.e. for <j) — ► 0, the diffusion tensor 
becomes diagonal and independent of r N , 

Dy(r^) =D 16ij (2.5) 

with 1 the unit tensor and 5ij the Kronecker symbol. However, in con- 
centrated suspensions, where hydrodynamic interactions no longer can be 
neglected, Dy(r^) becomes a function of the positions of all particles, in- 
volving therefore many-particle interactions. 

The diffusion tensor T)ij(r N ) is directly related to the experimental short 
time self-diffusion coefficient DA 6) 

D s (<f>) =< Q • T>u(r N ) ■ Q > eq (2.6) 

for any particle i. Here Q is a 3-dimensional unit vector. D s (<p) reduces in 
the dilute limit to the Stokes-Einstein diffusion coefficient Dq (cf.eq.(2.5)). 
D s (<f)) is a purely hydrodynamic quantity, which involves the calculation of 
the very complicated many-particle interactions, and has been the subject of 
research for many years, both theoreticailyO'~!IIl and experimentally' 1£H24 



By now the behavior of T){j(r N ) and D s (4>) for intermediate volume fractions 
up to <p ~ 0.45, is fairly well understood theoretically. For higher concen- 
trations to the best of our knowledge only semi-empirical results exist. 

However, in this paper we are particularly interested in the high volume 
fractions (f> > 0.40. Therefore we were forced to incorporate hydrodynamic 
interactions in an approximate (mean-field-like) fashion, using eq.(2.6). This 
approximation seems justified for high frequencies, where the particle distri- 
bution in the suspension is very close to the equilibrium particle distribution, 
so that the hydrodynamic interactions are described in first order by the hy- 
drodynamic interactions of the suspension at infinite frequency. The mean- 
field approximation will be done explicitly in the next section. For D s {4>) we 
use at the end a semi-empirical relation which is consistent with experiments 
and hard sphere computer simulations for all up to ~ 0.6oBrt. 

3 Large frequency viscosity for soft potentials 

We will proceed with the calculation of the asymptotic behavior of r]((f), uj) 
for large to by calculating the stress-stress auto correlation function p v (t; 4>) 
of eq.(2.2) with the microscopic stress tensor of eq.(2.3) and the iV-particle 
Smoluchowski operator of eq.(2.4) for a soft, but very steep potential. 

We restrict ourselves to pairwise additive potentials, i.e. $>(r N ) = S^ =1 l/(rjj), 
with V(rij) the two particle potential and = rj — r^, rjj = |rj,-|. We can 
then write eqs.(2.3) and (2.4), respectively as 

n y {r N ) = -\hiJ-^ 1 (3-D 
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and 

tyr" ;</>) = J2 (Sik^i-P(l-S ik )^^)--D ij (r N )-V j (3.2) 

i,j,k=l ^ * ' 

We determine the asymptotic behavior of 77(0,0;) from eq.(2.1) for a soft, 
but steep pair potential V\(r) = e(a/r) 1 , with e the pair interaction energy 
at r = a and I = I rJ^lnVJ(r)| the steepness of the potential. Since the 
hard sphere potential is the limit for I — > 00 of Ve(r), one can obtain the 
hard sphere result, by letting I —* 00 at the end of the calculation. This is 
discussed in the next section. 

In order to compute p„(t;<f)) for short time, we first expand /t>«(t;0) for 

Vj(r) in powers of t. Thus, we write first eq.(2.2) for p n (t\<f>) as 

00 

/^) = E-*X(^) ( 3 - 3 ) 

n=0 n - 

and then calculate p£(0) =< E^r^JP^; (f))H% y (r N ) > eq , using eq.(3.1) 
for the second E^L(r ) within the brackets, i.e.: 

^(0) = - 1) < ^ xy (r N W(r N ; 0)r 12 , x ^-^ > eq (3.4) 

for the soft potential V/(r). 

In this paper we will restrict ourselves to only the leading order in the 
steepness. Thus, we take into account, for each n, only the most divergent 
terms in /. This implies that, for short times, we can neglect all the contribu- 
tions of more than two particles to the equilibrium ensemble average < > e q, 
i.e., we only have to calculate the two-particle contributions, involving the 
equilibrium radial distribution function g eq (r;(J)). This, because, relative to 
the two particle contributions, the n-particle contributions are of the order 



/ ( n 2 \ so that they can be neglected in the limit I — > 00 



Thus, for Q(r N ; (p) in eq.(3.4) we can use eq.(3.2) with only i, j G {1, 2}, 
giving 



p n M) 



l -N{N-l) 



< E% y (r N ) [(Vi - • V r (r N ) ■ Vi 



r 12, 



> 



(:>(/ 



(3.5) 



Here we used the symmetry of fi(r ; cf>) in the particles 1 and 2, when applied 
on functions of ri2- We have introduced the relative diffusion coefficient of 



two interacting spheres D r (r N ) = 2(Du(r N ) — Di2(r Ar ))o In the dilute 
limit for just two particles, distant from all others, the diffusion tensors 



Dn(ri,r2) and Di2(ri, r2) are known 1111' nip c§ I . 

For concentrated suspensions we make a mean-field approximation. We 
replace ~D r (r N ) in eq.(3.5) by its mean value < D r (r N ) > eq , which reduces 
for high frequencies, i.e., for short times to twice the single particle short 
time self-diffusion coefficient D s (cp) as given in eq. (2.6)111. Thus we 
write in eq.(3.5), 



V r (r N ) =< D r (r JV ) > eq = 2D S {4>)1 



(3.6) 



consistently with eq.(2.6). 

Using this approximation and eq.(3.1) for i,j £ {1,2}, i.e. neglecting 
again all but two-particle contributions, we find straight forwardly: 



1 



p^) = -N(N-l)(2DM) n 



< r i2 , 



■u 



Vi 



dVi(r 12 ) 



■ Vi 



9Vi(r 12 ) 

ri2,x^ >eg(3.7) 



Since Vir^V^r^) = ri2ViVJ(ri2)(l+0(Z *)), we can shift, to leading order 
in the steepness, the differential operator Vi through ri2 in any product of 
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ri2 and Vi(r 12) or its derivatives. Thus, from eq.(3.7) we obtain: 

p n v ( ( p)= 1 -N(N-l)(2DM) n 
2 2 

• < Tj prV{{r) [(V! - (5V{{r)v) ■ V 1 ] n V/(r) > eq (3.8) 

with r = ri2, r = |ri 2 |, f = r/r , dVi(r 12) / dri 2 = V((r). 

Changing to spherical coordinates, using the definition of the equilib- 
rium radial distribution function g eq (r\<f>) for the soft potential Vi(r) and 
performing the angular integration we find from eqs.(3.3) and (3.7): 

o r 00 

p r) (t;<f>) = ^n 2 V / ^(rj^r^'W^M^-W^^'M (3.9) 
15 Jo 

where V r = d/dr. 

Thus we have expressed p-nit; <p) at large volume fractions cfi and to lead- 
ing order in the steepness I in terms of a one dimensional integral over r 
involving the high density equilibrium radial distribution function g e q(r; eft) 
for a potential Vj(r) of finite but large I and the effective short time self- 
diffusion coefficient D s (cp). In the next section we consider the hard sphere 
limit (I — > 00) of eq.(3.9). 

4 Hard sphere limit 

To evaluate eq.(3.9) we introduce the function 

y eq (r;0)=g eq (r;0y v ^ (4.1) 

As discussed in ref.[pj, y eq (r;(p) is, unlike g eq (r;4>), a smooth continuous 
function of r for a// r and Z. For hard spheres g^ q (r; <fi) = y^ q (r; <fi) for r > a 
and x( < / ) ) — fleij ( r = a l 4*) = 2/eo ( r = a \ 4>) 1S the P arr correlation function at 
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contact. Writing Vh s [r) = V^^r), using that 

e ~PV hs{ r)y L[r) = _}_ 5{r _ a) ( 4. 2) 

and using eq.(4.1), we can write for the stress-stress auto correlation p^t] 4>) 
of eq.(3.9) in the hard-sphere limit / — ► oo: 

p,(^) = -™ 2 ^) j o drr 4 e-^^VUr)e^ at VUr) (4.3) 

Here we have defined fi(f a = 2D s (<j))(V 2 r - /3V^,(r)V V), the radial part of the 
two-particle Smolochowski operator for a hard-sphere potential in relative 
coordinates. 

This expression can be calculated by a method, similar to that of Ci- 
chocki and Felderhoffl . The actual calculation is given in more detail in the 
Appendix, giving: 



Eq.(4.4) leads, with eq.(2.1), to our final result for hard spheres: 

niM ™ tiM) + A{<i>)±±L m (4.5) 

with the coefficient of the square root singularity A{4>) given by: 

^) = l^w(4))' <«> 

This result, based on the Green-Kubo formula for the frequency dependent 
viscosity r/(<^>, cj)® of a colloidal suspension consisting of hard spheres with 
hydrodynamic interactions included, is compared with experiments of van 
der Werff et alH in fig.l. Here we have used for x(0) the Carnahan-Starling 



approximation!^' till x(</0 = (1 — O.50)/(l — (f)) 3 for (j) < 0.5 and a one- 
pole approximation© x(0) = 1-2/(1 — (p/fim) for (j) > 0.5, with <p m = 0.63 
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the volume fraction at random close packing. For D s (4>)/Dq we have used 
Beenakker and Mazur's expression [0 for < 0.45 and D s (<p)/Do = 0.85(1— 
<t>/<f>m) for <p > 0.45H. 

Figure 1 clearly shows that in order to obtain agreement with experi- 
ment it is neccesary to include hydrodynamic interactions, i.e., to take into 
account the diffusion tensor Dy(r^), in the basic N-particle Smoluchowski 
equation. In a mean-field approximation this leads to the replacement of 
the Stokes-Einstein diffusion coefficient Do by the short time self-diffusion 



coefficient D s (cj)), a replacement also made by Brady£3. Eqs.(4.5) and (4.6) 
reduce to the exact expression (eq.(1.5)) obtained before by Cichocki and 
Felderhof for low densities!!] to 0{4> 2 ) (see Appendix). 

We will show in the next section (eqs.(5.1) and (5.3)) that the right hand 
side of eq.(4.4) is the leading term of the expansion in powers of I of p(t; <p) 
for a soft potential, i.e., for finite I and for frequencies w up to ~ 7^ °Dq^ ■ 

5 Soft potential 

We note that the large w-behavior of rj((p,uj) ~ 1/y/uJ is typical for hard 
spheres. For any soft, but steep potential rj(cp,uj) ~ l/u for uj — » 00. For 
example, the presence of a lubrication layer causes a change in the relative 
diffusion of two spheres at very short times, which leads to X/oj behavior at 



very high frequencies, as discussed by Cichocki and 

FelderhofH 

and Ralli- 



son and Hinchi53. To study the transition from the l/w-behavior (for any 
finite I) to the 1/y^-behavior (J = 00) we have calculated p v (t; 4>) of eq.(3.9) 
for finite /. For finite I eq.(3.9) can no longer be calculated by a method, 
similar to that of Cichocki and Felderhof §, described in the appendix, but 
involves the calculation of the complete eigenvalue problem of the radial part 
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of r = 2.D S (</>) (Vj? — PV/ (r)V r ) , the two-particle Smoluchowski operator in 
relative coordinates, for finite 1^3. The result for p v (t; <j)) in eq.(3.9) can be 
written 

2im 2 Vla 3 x( 
15/3 2 

with r(» = 2D s ((p)tl 2 /cr 2 . The function r(x) can be expanded for x <C 1 
as: 

r(sc) = 1 - x + ^x 2 + 0(x 3 ) (5.2) 

and for x > 1 as@^ 

r W = ^( 1 + ^ + (?)) < 5 - 3 > 

Thus to leading order in the steepness Z, we obtain from eqs.(5.1) and (5.3): 
2irn 2 Vl<j 3 X (4>) 1 , tj> D 

consistent with eq.(4.4) when I — > oo. For finite I, r](<f),u)) ~ l/y/uj for 
frequencies u up to ~ ;J^f, while for lar ger u, r]((f>,uj) behaves as 1/cJ, 
as is typical for soft potentials. It might be interesting to see whether 
such a transition in the asymptotic behavior of rj((f), u) can be observed in 
concentrated colloidal suspensions, where the interaction potential is steep. 
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Appendix 



CO 



p v (t; 0) = - 2 ™ 2 *ffiM / d r5(r - a)e^X» (A.4) 



Here we calculate the stress-stress autocorrelation function p v (t;4>) for 
hard spheres. We start with eq.(4.3) for p v (t;cf)) in the hard-sphere limit, 
i.e. I — > oo: 

^(t; 0) = ^vrn 2 y X (0) drrV^)^ (r)e ^ V^(r) (A.2) 
with Vh s (r) the hard-sphere potential and 

n h r s = 2D s ^){Vl - pvL(r)Vr) (A.3) 
Using that e^.p(-f3V hs (r))V/ ls (r) = -^S(r - a) gives 

) 

15/3 Jo 
With eq.(2.1), eq.(A.4) gives 

o /-oo /-oo 

= r^) - -TrnVxW J q dt J q drS(r - a)e^+^ V > s{r) 

2 f 00 1 

= ??oo + — 7mV 4 x(0) ^ dr^(r-(T) + ~ (A.5) 

We define 

f ^ ) = WT^) VL{r) (A - 6) 

and deduce the following differential equation for f(r,uj): 

(2D s V 2 r - 2D s pvl s (r)V r + iu) f(r, u) = V^r) (A.7) 

Due to the singular behavior of the hard-sphere potential at r = a, eq.(A.7) 
reduces to the boundary value problem: 

/ (2D a V* + io;)/(r,u;) = r>a 
\ 2D s pV r f(r,u) = -1 r = a 
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with the solution, bounded for r > a 

«'■"> = dn,^ (A ' 9) 

where a 2 = —iuj/2D s . Eq.(A.5) with eqs.(A.6) and (A. 9) gives the asymp- 
totic result for the frequency dependent viscosity as given in eqs.(4.4) and 
(4.5) 

n(M = + f^xW (j^f^ (A.10) 

where we have used the Peclet time Tp = a 2 /ADq and the Stokes-Einstein re- 
lation for Dq as given in eq.(1.3). Eq.(A.lO) with eq.(2.1) yields for p v (t; (j>): 

For low concentrations (<fi — > 0), x(</>) = 1 and D s (<p) = Dq and p^(t;^>) 
reduces to the result of Cichocki and FelderhofH for pr)(t; 4> — ► 0), for short 
times. These authors calculate pn(t] (f> — * 0) for hard spheres on the basis of 
eq.(2.2) restricted from the beginning to two hard-sphere particles only, but 
for all times t. One can show that the angular dependences in eq.(2.2) (i.e. 
in T, xy ) are irrelevant for short times. Thus, for t — > 0, both approaches are 
similar, leading to identical results for p n (t\4>) and rj(<p,u). 
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Figure captions 



Figure 1. Coefficient A{(ft) of the square root singularity in u, i.e., ~ A((ft)(l+ 
i)r]o/ \Jljtp, of rj((f>,uj) as a function of the volume fraction (ft. Experimen- 
tal points (•) of van der Werff et al, ref.Q; The vertical lines indicate the 
estimated errors in the experimental values of A((ft), while the horizontal 
lines indicate the effect of the 4% uncertainty in (ft (ref . |35|] ) . The dashed 
line represents the mode-coupling result eq.(1.2) and the solid line our re- 
sult given in eq.(4.6). Here we have used for xi'P) t ne Carnahan-Starling 
approximation EH x(0) = (1 — O.50)/(l — (ft) 3 for (ft < 0.5 and a one-pole 
approximation ¥3 x(^>) = 1-2/(1 — 4>/(ft m ) for (ft > 0.5, with <ft m = 0.63 
the volume fraction at random close packing. For D s {(ft)/DQ we have used 
Beenakker and Mazur's expression 10 for (ft < 0.45 and D s ((ft)/Do = 0.85(1 — 
(ft/4> m ) for 4> > 0.45c3. The dotted line represents eq.(4.6) with D s ((ft) = Dq, 
i.e., when hydrodynamic interactions are neglected. 
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